Spontaneous interlayer superfluidity in bilayer systems of cold polar molecules 
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Recent experimental progress in producing ultracold polar molecules with a net electric dipole 
moment opens up new possibilities to realize novel quantum phases governed by the long-range and 
anisotropic dipole-dipole interactions. In this work we predict the existence of experimentally ob- 
servable novel broken-symmetry states with spontaneous interlayer coherence in cold polar molecule 
bilayers. These exotic states, which are manifestations of collective bilayer quantum entanglement, 
appear due to strong repulsive interlayer interactions and exhibit properties of superfluids, ferro- 
magnets and excitonic condensates. 
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During the last decade we have observed the spec- 
tacular progress in the realization of various quantum 
phases using cold atoms. This progress has deepened our 
understanding of various phenomena such as BCS-BEC 
crossover of fermions pQ and superfluid-Mott insulator 
phase transition of bosons in an optical lattice [2 . How- 
ever, the variety of quantum phases that can be realized 
in cold atom systems is limited by the short-range nature 
of the interparticle interactions. Recent progress in pro- 
ducing and manipulating heteronuclear polar molecules 
[3] provides an opportunity to realize a plethora of novel 
quantum phases of matter governed by long-range in- 
teractions [4]. This interesting prospect is made possi- 
ble by the fact that polar molecules have large electric 
dipole moments associated with their rotational excita- 
tions, which lead to strong, long-range and anisotropic 
dipole-dipole interactions. The interactions between such 
polar molecules can be tuned using dc and ac electric 
fields pQ. Below, we concentrate on one intriguing aspect 
of fermionic polar molecule systems - the possibility to 
realize bilayer superfluidity (or equivalently bilayer XY 
ferromagnetism) with spontaneous interlayer coherence 
using interlayer molecular repulsion. 

Because of the fermionic nature of the molecules even 
for the spinless fermion case considered here, the inter- 
action has an exchange component in the layer (or pseu- 
dospin) index which drives the instability towards the bi- 
layer superfluid phase. The dominant contribution to the 
exchange energy comes from the short distances (large 
momenta) where the interlayer dipolar interaction is re- 
pulsive. Thus, the novel collective bilayer state we pre- 
dict arises from a repulsive interaction in sharp contrast 
to all other superfluid quantum phases discussed in cold 
atomic fermions where inter-particle attraction leads to 
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superfluidity. In the symmetry broken interlayer coher- 
ent phase, the particle number in each layer becomes in- 
determinate in spite of the interlayer single-particle tun- 
neling amplitude being almost zero. Such a state is very 
analogous to an excitonic superfluid in which excitons 
formed a quasiparticle in one layer "binding" to a quasi- 
hole in the other layer condense into a phase coherent 
state. Up to this date the clearest evidence for the real- 
ization of this type of exciton superfluid state has been 
observed in semiconductor bilayers in the Quantum Hall 
(QH) regime, in which the layers are immersed in very 
high magnetic fields 0. The unavoidable presence of 
disorder in solid state systems as well as the nature of 
the measurement involving finite interlayer tunneling [8 
cause complications in the pristine realization of the in- 
terlayer superfluid phase in QH bilayers. As a result, 
vortices and the Berezinskii-Kosterlitz-Thouless (BKT) 
transition, unambiguous signatures of the interlayer co- 
herent state, have not been observed yet in solid state QH 
systems. The predicted polar molecule bilayer superfluid 
phase should be more striking because of the lack of disor- 
der in cold atom systems, the tunability of the interaction 
strength and the availability of experimental techniques 
allowing imaging of vortices [9 . Moreover, by adding 
an optical lattice potential in the xy-plane it is possible 
to modify the single-particle dispersion of the particles 
and model various condensed matter systems. For ex- 
ample, the bilayer cold-polar-molecule system with hon- 
eycomb lattice potential will mimic the bilayer graphene 
system. Thus, the realization of this interesting phase in 
cold atom systems is of great importance for understand- 
ing the instabilities driving the bilayer superfluidity and 
in general the physics of exciton condensation. 

Theoretical model Our starting point is the Hamilto- 
nian for fermionic polar molecules tightly confined along 
the z-direction by the laser field as shown in Fig. [T}d. 
We consider two clouds of polar molecules separated by 
a distance l z much larger than the confinement length 
w z of the molecules within each layer. When the con- 
finement length w z is much larger than the size of the 
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FIG. 1: (Color online) (a) Heteronuclear polar molecules AB 
with AB being 40 K 87 Rb, 7 Li 40 K, 6 Li 133 Cs. (b) Bilayer sys- 
tem of cold polar molecules in the presence of a circularly 
polarized ac electromagnetic field propagating along the z- 
direction. Schematic picture of different phases in the bilayer s 
system: normal (c) and pseudospin ferromagnetic (d). 



polar molecules, the rotational motion of the molecules 
is three-dimensional and is described by a 3D rigid ro- 
tor Hamiltonian. Polar molecules have permanent elec- 
tric dipole moment d, which couples to internal rota- 
tional degrees of freedom. The dipole moment leads to 
long-range interlayer and intralayer dipole-dipole interac- 
tions. The Hamiltonian of the polar molecules H reads 
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FIG. 2: (Color online) (a) Intralayer Born-Oppenheimer po- 
tential for polar molecules. The dashed(blue), solid(red) 
and dot-dashed (violet) lines correspond to Qr/A = 1/4, 
£}r/A = 1/8 and Qr/A = 1/20, respectively, (b) Interlayer 
Born-Oppenheimer potential for Qr/A = 1/8. The dash- 
dot (blue), solid(red) and dashed(brown) lines correspond to 
Iz/pA = 3,2, 1.5, Iz/pA — 1.5, respectively. For typical inter- 
particle distances considered here the interlayer interaction is 
repulsive. Here Qr/A = 1/8 and l z /pA = 3. 



of the microwave field. In this case, the effective inter- 
action between polar molecules V e ^(r) is given by the 
dressed Born-Oppenheimer potential adiabatically con- 
nected to the state see Fig.l of the supplemen- 
tary material. At large distances the dipolar interaction 
can be obtained perturbatively by first calculating the 
effective dipole moment (+|d|+) = d e ^ (coscjt, sincjt, 0) 
with d e ff = — y/2a+a-d t . The time- averaged interaction 
between dipoles in layers A and A' takes the form 
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where p = (p x ,P y ) is the center-of-mass momentum of 
a molecule with mass m, is the distance between 
two molecules, B is the effective rotational energy and 
J = (J x , J y ,J z ) is the angular momentum operator. The 
rotational eigenstates are \J,Mj) with J and Mj denot- 
ing the total internal angular momentum and its projec- 
tion on the quantization axis, respectively. 

The orientation of the dipole moments can be con- 
trolled with dc and ac electric fields Fig.[T]3, pQ . The tran- 
sition dipole moment between the states with J = and 
J = lisd t = | (0,0|d|l,Mj) | = d/y/3 with Mj = 0,±1. 
A circularly polarized ac electric field E ac (t) propagating 
along z-direction, see Fig. []J>, drives transitions between 
the rotational states |0, 0) and |1, 1) with Rabi frequency 
Qr = dtE ac . If the frequency of the field uo is close to the 
transition frequency ojq = 2B between the states |0,0) 
and |1, 1) (i.e. the detuning A = oj — ujq <C ojq), the lead- 
ing effect of the electric field is to mix these two states. 
Within the rotating wave approximation, the dressed 
states are given by |±) = a±|0,0) ± a T e~ luJt |1, 1), 
where a + = -T/y/T* + Sl% a_ = Vt R / ^T 2 + W R and 
2T = A + s/A 2 + 4ft| [DHl. Polar molecules can be pre- 
pared in the internal state |-f) i by an adiabatic switching particles is given by the dressed Born-Oppenheimer po- 



where p = (x, y) is the 2D coordinate, and z\\> = l z for 
A ^ A' and zero otherwise. At short distances, when the 
dipolar interaction energy is comparable with the detun- 
ing, the above perturbative treatment breaks down. In 
order to find the Born-Oppenheimer potential at short 
distances p < Pa = (d 2 /A)z, it is necessary to account 
for all couplings between different angular momentum 
channels within the J = 0, 1 manifold pQ. The exact 
Born-Oppenheimer potentials are shown in Fig. [2] One 
can notice that the effective intralayer dipole-dipole in- 
teraction between polar molecules prepared in the state 
|+) becomes repulsive at p ~ pa clue to the presence 
of avoided crossings with other field-dressed levels [4]. 
For 6 Li 133 Cs molecules, typical parameters are d « 6. 3D, 
B w 6GHz, and A & 10MHz yielding the length scale 
Pa ~ 50nm that is much larger than the characteristic 
scale of dipole-dipole interactions p# = (d 2 /B) 1 ^ 3 ~ lnm, 
which sets the short-range cutoff. Thus, for typical den- 
sities of polar molecules considered here, no ~ 10 T cm -2 , 
the ac electric field shields the molecules from short-range 
inelastic collisions and prevents from the collapse of the 
system P3HJ[10], see supplementary material. 

Henceforth, we consider the dilute gases of polar 
molecules, where the interparticle distance is larger than 
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tentials shown in Fig. [2j In order to avoid unwanted in- 
elastic collisions leading to the decay of the molecules in 
an s-wave channel, we assume the molecules to be spin- 
polarized. In this limit, the effective second-quantized 
Hamiltonian of the bilayer system takes the form 
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where cj^x an d c\ x are the fermion creation and anni- 
hilation operators for a molecule with momentum k in 
layer A. The strength of the dipolar interactions can 
be characterized by the dimensionless parameter r s = 
dg ff m 010/271-. As r s is increased, the bilayer system be- 
comes susceptible to various instabilities driven by the 
dipolar interactions. Here we concentrate on the insta- 
bilities induced by the interlayer interactions and neglect 
the instabilities induced by the intralayer interactions, 
see supplementary material for the justification of such 
approximations. 

To understand the interlayer instability, it is conve- 
nient to draw an analogy with ferromagnetism and intro- 
duce pseudospin-1/2 operators 2m^ = c\ x cr w Ci\> . The 
spinors \\) and |^) represent the states in which molecules 
are in layer 1 or 2, respectively. When l z ^> fc^ 1 , with 
kp = ^Attuq being the Fermi momentum, the molecules 
in different layers are uncorrelated, and the many-body 
state of the system is given by 
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For equal densities in the layers ni = n 2 = no, the total 
magnetization M is zero, M = (^n| J^i^i I^n) = 0, 
similar to paramagnets. The normal state |^n) mini- 
mizes the kinetic energy at the expense of the potential 
energy, which is at its maximum. When the interlayer 
distance becomes smaller than fc^ 1 , the potential energy 
becomes large and at some point starts to dominate over 
the kinetic energy. In this case, the system favors the 
state in which fermions in different layers are correlated 
in a way that minimizes the interaction energy, i.e. the 
molecule in layer 1 is coupled to a "hole" in layer 2. At 
the mean-field level such correlations are captured by an 
order parameter A12 ex (cj^c/^) 7^ 0. Since the prod- 
uct wavefunction does not have such entanglement 
between the layers, the bilayer system should undergo 
a quantum phase transition as a function of the dis- 
tance l z or the strength of the dipole moment <i e ff- The 
many-body wavefunction minimizing the interaction en- 
ergy takes the form 
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In this entangled state the state of the molecule is given 
by the coherent superposition of the amplitudes in dif- 
ferent layers. Thus, even in the absence of tunneling, 
the molecule layer index becomes uncertain. Using the 
spin analogy, the state |^fm) has non-zero magnetiza- 
tion M =(^FME i rn i |^FM) ^ with M lying in the 
xy-plane. Similar to superconductors, this ferromagnetic 
state |^fm) spontaneously breaks U(l) symmetry and 
develops interlayer coherence. In this state the phase 
difference between different layers tp is well defined and 
the number of molecules in each layer fluctuates satisfy- 
ing the uncertainty relations Am z Acp > 1/2 [11] despite 
the absence of interlayer tunneling in Eq. ([3|. 

The phase diagram between the two competing states 
- normal |^n) and pseudospin ferromagnetic |^fm) - can 
be obtained using variational mean-field calculation [12]. 
The total energy per area A of the bilayer system in the 
normal phase is 
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where the dimensionless function Z(a\ which describes 
intralayer interaction, is defined as 
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with f(x) — arccos(x) — x\J\ — x 2 . The momentum de- 
pendence of the intralayer interaction potential V^\q) 
is shown in Fig. 1 (b) of the supplementary material. 
One can notice that the contribution of the interlayer 
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(4) Hartree term is zero here because V^(q) = d^irqe' 



goes to zero as q — » 0. To calculate the energy of the 
system in the interlayer coherent state |^fm) we first 
introduce the order parameter 



Ai 2 (fe) 
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which takes into account interlayer correlations. The or- 
der parameter is obtained by numerically solving the self- 
consistent equation above, subject to the total particle 
number conservation constraint. Because of the dipolar 
nature of the interaction Ai 2 (k) has momentum disper- 
sion. The dependence of Ai 2 (0) on d e ^ and l z obtained 
self-consistently is shown in Fig. [3] (a). At the mean-field 
level, the Hamiltonian ([3| can be diagonalized using a 
Bogoliubov transformation yielding the many-body vari- 
ational wavefunction ([5|. For sufficiently large interac- 
tions the lowest energy state of the system corresponds 
to the pseudospin ferromagnetic state fully polarized in 
the XT/-plane [12 , see also supplementary material. The 
total energy of the system per area in the interlayer co- 
herent state is given by 
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FIG. 3: (Color online) The dependence of Ai2(0) and AE on 
l z and deff is shown in (a) and (b), respectively. Here we used 
pA — 50nm, no — 10 7 cm -2 and m — 139ra p with m p being 
the proton mass. The dependence of AE on no and d e ff at 
fixed l z — 0.3/xm and on no and l z at fixed d e s = 1.5 D is 
plotted in (c) and (d), respectively. 



where So = ^im^jm and the function F(a) describing 
the interlayer exchange contribution reads 



F(a) = 32V2tt J x 2 dx ^arccos(x) — xy/l — x 2 



(10) 



The energy difference between normal and ferromagnetic 
phase AE/A — (Eym. — F^)IA determines the mean- field 
phase diagram for the bilayers of polar molecules shown 
in Fig. [3} 

The long wavelength Hamiltonian describing the phase 
fluctuations is 



H= \ I d 2 r Ps \V^\ 2 , 
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where p s is the "spin stiffness", which is the result of 
the loss of interaction energy due to the spatial vari- 
ations of order parameter phase ip. The effective XY 
model defined by Eq. (11) undergoes BKT transition as- 



sociated with unbinding of vortex pairs at the temper- 
ature Tbkt ~ ttPs/2. The "spin stiffness" can be cal- 
culated within linear response theory, which yields the 
result p s = no /2m, similar to the one in superfluids. 
Thus, the BKT transition in the bilayer system of po- 
lar molecules occurs at the temperature Tbkt ~ £p/8 
with sp — kp /2m, that for no ~ 10 7 cm -2 corresponds 
to a temperature that should be accessible in the near fu- 
ture [13 . In the pseudospin ferromagnetic phase vortices 
are correlated in different layers, and the BKT transition 
in the bilayer system can be detected by imaging vortices 
using matter wave heterodyning techniques [9]. 

In summary, we predict an unusual broken-symmetry 
phase with spontaneous interlayer coherence in a bilayer 
system of cold polar molecules. Our main findings, sum- 
marized in the phase diagram shown in Fig. |3j indicate 
that the experimental observation of such phase requires 
low densities of cold polar molecules no ~ 10 T cm -2 , re- 
alistic dipole moments <i e ff ~ ID, and reasonably low 
temperatures T ~ InK. Given that for these parameters 
the inelastic decay rate is small [4 , the detection of this 
exotic many-body state should be within the experimen- 
tal reach in the near future. The novel phase we predict is 
an interlayer entangled state, arising from the repulsive 
part of the dipolar interaction and exhibiting superflu- 
idity (or, equivalent ly, XY pseudospin ferromagnetism) 
between the layers rather than within the layers. 
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I. SUPPLEMENTARY INFORMATION 

In this document we elaborate on the analysis of the 
various instabilities in the bilayer system of cold polar 
molecules. We begin with the discussion of the role of 
the intralayer interactions in stabilizing the cold-polar- 
molecule Fermi gas. 



A. Stability analysis of the single layer 

The Hamiltonian for a single layer 2D Fermi gas of cold 
polar molecules reads 



!hb 



k+q^k' 
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(12) 

The intralayer interaction is given by the dressed Born- 
Oppenheimer potential Vu(p) [HE], which has a short- 
distance repulsive core and an attractive tail at large dis- 
tances as shown in Fig.[5^i of the main text. The crossover 
between these two regimes occurs at the length scale pa, 
which is determined by the detuning A, see Fig.|4j Thus, 
by changing the detuning frequency one can control the 
strength of the intralayer interaction. 

The energy of the single layer in the Hartree-Fock ap- 
proximation [3] is given by 
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where nk = Q(k F — k) with Q(x) being the Heaviside step 
function. The dimensionless function I(q, x) is given by 
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and the momentum dependence of the intralayer poten- 
tial A) (0) - V^ X \q) is shown in Fig. [5)0. Then, the 
total energy of the system becomes 
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with r s = cZg ff my / no/27r. The dimensionless function 
Z(a) is defined as 
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FIG. 4: Two-body energy levels as a function of the distance 
for ft*/ A = 1/8. 



where f(x) = arccos(x) — x\J\ — x 2 . Since the inte- 
gral ( [l6| ) involves both the contribution from the attrac- 
tive tail and repulsive core of the potential V^ 1 \p) (see 
Fig [5]), the intralayer interaction energy depends on the 
parameter k F p/\. For typical densities considered here 
no ~ 10 T cm -2 and pA ~ 50nm (i.e. k F pA ~ 0.06), the 

compressibility of the system ft -1 = ^ Q ^ 2 is positive, 
see Fig [5^, and thus, the ac electric field allows one to 
stabilize the Fermi gas and prevent phase separation. 



B. Intralayer superfluid instability. 

We now discuss the possibility of the superfluid insta- 
bility of the Fermi surface in a single layer of cold polar 
molecules. It was shown in Ref. [4 that the intralayer in- 
teraction can lead to p-wave pairing. The weak coupling 
BCS theory developed in Ref. [4 is valid for r s <C 1, 
and the evaluation of the superfluid T c for typical pa- 
rameters considered here (r 3 < 1) requires developing a 
strong coupling theory of superfluidity, which is beyond 
the scope of the present paper. However, we note that 
by tuning the ratio Qr/A the depth of the attractive po- 
tential can be decreased as shown in Fig. |5^i indicating 
that it is possible to substantially suppress the intralayer 
superfluid T c or eliminate such instability altogether. To 
demonstrate this, we calculate the intralayer interaction 
in the p-wave channel V^ =1 ) 

POO 

V (i) = V (l) (k F ,k F ) = 2tt / pdpjf{k F p)V {11 \p), 

Jo 

(17) 

which can be controlled by changing Qr/A. Here 
Ji(k F p) is the Bessel function of order I. Due to the 
short range repulsion the potential V^ =1 ) changes sign 
as a function of k F p&, see Fig. [5]i. For the relevant pa- 
rameters considered here Qr/ A = 1/8 and k F p& ~ 0.06, 
the potential V^ l=1 ^ is positive, and the analysis of the 
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FIG. 5: (a) Intralayer Born-Oppenheimer potential for 
polar molecules. The dashed (blue), solid (red) and dot- 
dashed (violet) lines correspond to different Rabi frequencies: 
Qr/A = 1/4, Qr/A = 1/8 and tt R /A = 1/20, respectively, 
(b) Fourier transform of the intralayer Born-Oppenheimer po- 
tential for VLr/ A = 1/8. (c) The compressibility of the single- 
layer cold polar Fermi gas as a function of density no- Here 
Ko = m/ir. (d) The dependence of the p-wave potential V^ l=1 ^ 
on the dimensionless parameter k F pA- 



spontaneous bilayer coherence neglecting the intralayer 
instabilities in the particle-particle channel is justified. 
We note here that there is no inconsistency with the con- 
clusion of Cooper and Shlyapnikov [4 since for the pa- 
rameters used in Ref. [4 (Qr/A = 1/4 and pA = 30nm) 
the potential is V^ l=1 ^ is indeed negative leading to the 
p-wave superfluidity. 



C. Effect of intralayer interactions on the 
pseudo-spin instabilities 

In this section we provide details on the calculation of 
the energy of the bilayer system of cold polar molecules 
within variational mean-field theory. We begin by an- 
alyzing the total energy of the system per area in the 
normal phase described by the wavefunction |^n) (see 
Eq. (4) of the main text): 
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The energy E^ essentially follows from Eq.(15) multi- 
plied by the number of layers. 

We now discuss the energy of the bilayer system in the 
superfluid state. Assuming that the dominant contribu- 
tion to the energy comes from interlayer interactions, we 
introduce the mean field Ai2(fc) (see Eq. (8) of the main 
text) and neglect for now the intralayer interactions. At 



the mean-field level, the Hamiltonian 
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can be diagonalized via canonical transformation yield- 
ing the energy spectrum e±(k) = (J^ — /i) ± |Ai2(fc)|. 
Assuming that all molecules occupy the lowest energy 
band S- (fc), the wavefunction corresponding to such state 
I^fm) is given by Eq. (5) of the main text. We calcu- 
late the order parameter self-consistently by enforcing 
the above constraint. Using the variational wavefunction 
I^fm), we then calculate the total energy of the system. 
The contribution of the intralayer interaction reads 
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and the total energy in the interlayer-coherent superfluid 
phase is given by 
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^ Jo 



1 - r s F(V8k F l z ) 

v£ 1 \0)-V£ 1 >(y/8k FPA x) f(x) 



m 



1 - r s F(Vsk F l z 



Z(V8k F PA ) 
k F p A - 



Here the function F(a) is defined as 

F{a) = 32V2nJ x 2 dx (^arccos(x) —x\J l-x 2 ^je~ ax 

[a(37r-4a + 157r/ (a))-37r(12 + a 2 )/i(a) 

(22) 



32V27T 



3a 4 

+3tt(2 + a 2 )Li(a) - 157raL 2 (a)] 



with I n (x) and L n (x) being the modified Bessel and 
Struve functions, respectively. As follows from Eq.(21), 



the magnitude of the intralayer and interlayer interac- 
tions depend on different physical parameters k F l z and 
k F pA- Thus, one can control their relative contributions 
to the total energy so that the isospin ferromagnetic state 
is favorable. 

We conclude this section by calculating the energy of 
the system in the Ising ferromagnetic state described by 
the wavefunction \I>im = Ylk<V2k F c ik 1^)- Such state has 
the pseudospin magnetization M aligned along the z- 
axis, which corresponds to all molecules being in a single 
layer. The intralayer interaction energy in this state can 



7 



be written as 



(IM) 



A 



k,q 



(23) 



For typical parameters considered here, the Ising pseu- 
dospin ferromagnetic state |^im) has higher energy than 
the state |\&fm)) and thus, can be ignored. 



d 2 i 



(2tt)2 L 



v£ x \0)-V^ x \ q )]l( q ,V2k F ) 



= -n 2 [ xte\v& 1 \0)-v£\>ftk FP ±xif{x). 
^ Jo 1 J 

The total energy of the system in the Ising ferromagnetic 
state reads 



£ t (IM) _W 32 



A 



m 



o 32 of 7 
-H n / xax 

n Jo 

z(Vsk FPA ) 



v£ 1 \0)-v£ 1 \V8k FPA xj\f(x) 



(24) 
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